#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Sat Nov 11 23:35:32 2023

@author: liqingsimac
"""
import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt

x=np.array([0.9,1.1,1.8,2.3,3.0,3.3,4.0])
p=np.array([0.37,0.31,0.44,0.60,0.67,0.81,0.79])

X = sm.add_constant(x)
md = sm.OLS(p,X).fit()
beta01A=md.params
print(beta01A)

beta01B=np.linalg.inv(X.T@X)@X.T@p
print(beta01B)

y = np.log(p/(1-p))
beta01C = np.linalg.pinv(X)@y
print(beta01C)

beta01D=np.linalg.inv(X.T@X)@X.T@y
print(beta01D)

fig=plt.figure()
ax=fig.add_subplot(111)
ax.plot(x,p,'bo')
ax.set_xlabel('x')
ax.set_ylabel('p')

xx = np.linspace(start=-2,stop=7,num=61,endpoint=True)
pp1 = beta01A[0] + beta01A[1]*xx
pp2 = 1/(1+np.exp(-beta01C[0]-beta01C[1]*xx))
ax.plot(xx,pp1,'r--',label='linear regression')
ax.plot(xx,pp2,'g--',label='logistic regression')
ax.legend()

# 输出结果：
# [0.17811832 0.16726657]
# [0.17811832 0.16726657]
# [-1.40066042  0.73900674]
# [-1.40066042  0.73900674]

